Magnus Force in Discrete and Continuous Two-Dimensional Superfluids 
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Motion of vortices in two-dimensional superfluids in the classical limit is studied by solving the 
Gross-Pitaevskii equation numerically on a uniform lattice. We find that, in the presence of a 
superflow directed along one of the main lattice periods, vortices move with the superflow on fine 
lattices but perpendicular to it on coarse ones. We interpret this result as a transition from the 
full Magnus force in the Galilean-invariant limit to vanishing effective Magnus force in a discrete 
system, in agreement with the existing experiments on vortex motion in Josephson junction arrays. 
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I. INTRODUCTION 

There is a long-standing interest in how the Magnus 
force, acting on vortices in superfluids and supercon- 
ductors, changes as the system moves away from the 
Galilean-invariant (GI) limit. Indeed, experiments in- 
dicate that the effective Magnus force is very small both 
in conventional bulk superconductors — except for very 
clean onesm — and in "discrete" superconductors, formed 
by Josephson-junction arrays (JJA)i^ In the first in- 
stance, there is a convincing explanation for this small- 
ness, based on the spectral flow of fermions at the vortex 
corei ?' 4 '^i7 The spectral flow creates an additional force 
on the vortex that reduces the total, effective Magnus 
force nearly to zero. However, the second case has re- 
mained something of a mystery. 

Various explanations of the smallness of the Magnus 
force in JJA have been reviewed in Ref. |8|. One proposal 
is that in this case the Magnus force is proportional not 
to the total density of electrons, but only to the "off- 
set charges" , given by the deviation of the system from 
electrical neutrality. 9 Another proposal is that the effec- 
tive Magnus force vanishes exactly as a consequence of 
the particle- hole symmetry. 10 However, Volovik 8 has ar- 
gued that the particle-hole symmetry in these systems 
is not exact, and as a result the effective Magnus force 
is nonzero, although small. Finally, we mention that 
when the Josephson barrier is metallic, cancellation of 
the Magnus force can be explained^ by a spectral flow 
mechanism similar to that in bulk superconductors. This, 
however, does not work when the barrier is insulating. 

There is a common theme to the above proposals: they 
all make use of specific properties of the electronic spectra 
or, alternatively, of the particle-hole symmetry already at 
the level of an effective description — in terms of phases 
and charges of superconducting islands. 10 As we will see 
below, that symmetry, present in the simplest model of 
JJA, results in fact from neglecting the coupling between 
phase gradients and density fluctuations. The question 
then is whether this assumption indeed applies in the dis- 
crete limit, or a nonzero Magnus force persists no matter 
how discrete the system becomes. 

In this paper, we report results of a numerical study 
of the Magnus force. These results have been obtained 



by numerical solution of the classical Gross-Pitaevskii 
(GP) equation in two dimensions (2D). Since the clas- 
sical approximation neglects the commutator of the Bose 
fields ^P, W in comparison with the average density 
n ave = it requires that the number of particles per 

site be large enough. More precisely, the classical limit in 
JJA is reached when the Josephson energy is much larger 
than the charging energy^ An equivalent condition is 
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where a is the lattice spacing (we assume a square or a 
nearly square lattice), and £ is the "healing" length, de- 
fined below. The right-hand side of (JIJ is a dimensionless 
measure of the interaction strength. 

To explore the role of discreteness, we solve the time- 
dependent GP equation on uniform spatial lattices with 
different values of the lattice spacing a. The relevant 
length scale to which a can be compared is the "healing" 
length £. When a < £, we reach the nearly GI limit, 
in which the GP equation describes a quasi-continuous 
neutral superfluid. Vortices have a core of size £, which 
is resolved by the lattice. In the opposite limit, a ^> £, 
vortices have no core, in the sense that there is no sig- 
nificant depletion of density anywhere. The lattice sites 
can then be thought of as corresponding to individual 
islands, each of which is characterized by a value of the 
phase variable — a model of a JJA. 

More precisely, the Lagrangian of our model in rescaled 
variables, for the case of a square lattice, is 



£ = XM { idt ~ o ^ i 2 + const ) ^ 



1 



a 



(2) 

where the first sum is over all lattice sites, and the 
second — over all nearest neighbors. Note that if we 
write ip = y/nexp(i9) and neglect fluctuation of density 
Sn = n — n a ve in the second (gradient) sum in Eq. J5J, 
the classical equations of motion become invariant un- 
der the transformation 9 — > —6, Sn — > — Sn. This is the 
particle-hole symmetry that was used in Ref. to argue 
the absence of Magnus force in JJA. Here, we study the 
complete Lagrangian (J2J) , including the coupling between 
Sn and the phase gradients. 
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On a coarse lattice, and in the presence of superflow, 
the rotational invariance is broken so strongly that defi- 
nition of the Magnus force becomes a non-trivial matter. 
Indeed, we have found that in general vortices and anti- 
vortices do not even move symmetrically with respect to 
the superflow. However, when the superflow is along one 
of the lattice's main periods, they do, and we concentrate 
on this case in what follows. 

We have found that while on fine lattices vortices move 
with the superflow, as expected in the GI limit, on coarse 
lattices they move perpendicular to it. We interpret this 
as vanishing of the effective Magnus force in the discrete 
limit, in agreement with the experiments on JJA. 2 

The paper is organized as follows. We describe details 
of the numerical procedure in Sect. |n| In Sect. II I II we 
describe a phenomenological model that we use to in- 
terpret our numerical results. This model allows us to 
convert measurements of the longitudinal and transverse 
velocities of a vortex into values of the effective Magnus 
and drag force coefficients. Numerical results are pre- 
sented in Sect. I1VI Sect. ]V\ is a conclusion. 

II. NUMERICAL METHOD 
A. Dimensionless variables 

The continuum GP equation has the form 

i™ = -±VH + 9m i % (3) 

where ^ is a complex scalar, the order parameter of the 
superfluid, m is the mass of a fluid particle (in units 
where % = 1), and g is the interaction constant. For 
computational purposes, it is convenient to scale out the 
parameters by expressing length, time, and the order pa- 
rameter in their "natural" units. A natural unit of length 
is the "healing" length 

£= (2m 5 n avc )- 1 / 2 , (4) 

where n avc = ^vptvp) is the average density of the fluid. 
Then, new, tilded variables are defined by the following 
relations: 

x = x£, (5) 
t = t^—, (6) 

fJTl&VC 

* = ^^/n^, ■ (7) 
In the new variables, the GP equation simplifies to 

i| = -VV+l^. (8) 
ot 

Note that by virtue of Q the rescaled average density is 
always equal to one: 

navo = (Vfy) = L (9) 
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FIG. 1: (Color online) A representative initial state with two 
vortices of topological charge 1 (+) and —1 (— ) before super- 
flow is turned on. 
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FIG. 2: (Color online) A representative initial state after su- 
perflow is turned on. 

B. Computational scheme 

To study the motion of vortices, we discretize Eq. (JSJ) 
on a uniform spatial lattice and solve it as an initial value 
problem, i.e., knowing the state of the system at some ini- 
tial moment, we calculate the state at later times. This 
requires imposing suitable initial and boundary condi- 
tions (see below). The lattice in general has different 
lattice spacings in the x and y directions. We used an 
operator-splitting algorithm with separate updates for 
the Laplacian and potential terms in JSJ. Updates cor- 
responding to the Laplacian were done using the Crank- 
Nicholson scheme, which is unconditionally stable. The 
complete algorithm is unitary and second-order accurate 
in space and time. 

C. Boundary and initial conditions 

To avoid effects of the boundaries on the motion of the 
vortices, we use periodic boundary conditions in both di- 
rections. The initial states for the runs are created in the 
following way. We begin with the following field, con- 
taining a vortex and an anti-vortex (the presence of an 
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anti-vortex is necessary to satisfy the boundary condi- 
tions): 

««>- n fe4 n ( 4^fr- »») 

z + e[z+] 1 +l *_e[z_] 1 1 

where Z + and Z_ are the desired (complex) positions of 
the vortex and anti- vortex, and [Z±] denotes the set of 
positions including Z± and a few mirror images with re- 
spect to the boundary. Then, evolving the system in 
the imaginary time, we cool the system down. Posi- 
tions of the vortices during the cooling do not change, 
so we can place the vortices in convenient locations. To 
minimize effects of the vortex- anti- vortex interaction, we 
place them half of the total lattice length apart. After 
that, we turn on a superflow and begin evolution in real 
time. A representative initial state, before and after the 
superflow was turned on, can be seen in Figs. ^ El 



u 




FIG. 3: (Color online) Diagram of forces acting on an 
(anti)vortex in an isotropic superfluid: FLorentz perpendicu- 
lar to the superflow velocity U, FMagnus perpendicular to the 
vortex velocity V, and Fdmg opposite to V . In dynamic equi- 
librium the sum of all forces vanishes, and the vortex moves 
at an angle 6 V = 7r/2 — Quail with respect to the superflow. 



D. Velocity measurements 

The aim of the simulations is to observe the motion 
of vortices in the presence of superflow. The order pa- 
rameter of the superfluid, being a complex scalar, can be 
written as 

tp = Vrie ie , (If) 

where n is the local density of the superfluid. Then, 2V# 
is the local superfluid velocity. The average velocity of 
the superflow U is calculated as an average over the entire 
lattice. The velocity of a vortex V, on the other hand, is 
obtained from direct tracking of the vortex position dur- 
ing the simulation. In a GI system, we expect vortices to 
move with the flow: V = U. If the invariance is broken, 
they may behave differently. To see the actual behav- 
ior, we break GI by solving the problem on increasingly 
coarser lattices. 



III. PHENOMENOLOGICAL MODEL 

While the vortex velocity can be measured directly in 
our simulations, converting these measurements into a 
value of the Magnus force requires a model of forces act- 
ing on the vortex. A fairly conventional model is available 
for an isotropic fluid (which we expect to apply also in the 
GI limit on a lattice), but on coarse lattices modifications 
are needed. In this section, we review the conventional 
model, and then describe new effects introduced by the 
lattice. 



A. Magnus force in isotropic fluid 

The conventional (phenomenological) model includes 
three forces acting on a vortex, see Fig. [31 First, there 



is a drag force Fd ra g = —^V, directed against the vortex 
velocity. It accounts for longitudinal momentum transfer 
from the vortex to the lattice and to excitations (quasi- 
particles). The latter channel is dissipative, i.e., the work 
done by the drag force goes into excitation of the quasi- 
particle subsystem. Close to the GI limit, we expect the 
classical momentum transfer to quasiparticles to be inef- 
fective, and hence the drag to be small. 

Second, there is an effective Magnus force FMagnus = 
—<t v z x V , perpendicular to the vortex velocity. 

Finally, there is a force, perpendicular to the superflow 
velocity and accounting for the work done by vortices as 
they unwind the superflow. We refer to it as the Lorentz 
force, FL oren t z = a u z x U. The coefficients i], a u , and a u 
refer to unit inertial mass. We will never need to discuss 
the actual value of the vortex inertial mass in this paper. 

So, the equation of motion for the vortex is 

dV _ - 

^ — Fdrag H~ -^Magnus ~T" FL oren i z 

= —r)V — a v z x V + a u z x U . (12) 

In complex notation, where we identify the ^-direction 
with the real axis and the y-direction with the imaginary 
axis, the equation becomes 

dV , x 

— = -aV + (3U, (13) 

where a = rj + ia v and (3 = ia u . 
The solution is easily found to be 

V = -U + C exp(-at) , (14) 
a 

where C is an integration constant. The exponential term 
is a transient that rapidly decays and turns out to be too 
small to be observed even at small times. At large times, 
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it drops out altogether. Then, the solution becomes a 
motion with a constant velocity at an angle 9 V , given by 







tan#„ = arg — = rj/a v , 
a 



(15) 



with respect to the supercurrent. The angle 9 V is related 
to the Hall angle 0Haii frequently used in the literature 
by 9 V = 7r/2 - 9 Hail- 

The conclusion that vortices move in straight lines is 
well born out numerically. Notice that the steady velocity 
V = (j3/a)U depends only on the ratios a v fa u and r//<r u , 
characterizing the Magnus and drag forces. Measuring 
two components of the steady velocity, we obtain two 
equations for these two ratios, which can be solved with 
the result 



VxUx + VyUy 



VyU X 



V 2 



V X Uy 



v' + v 2 



(16) 
(17) 
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FIG. 4: (Color online) Longitudinal and transverse velocities 
of the vortex in units of the superflow velocity for different 
lattice spacings. The length of the lattice is 600; the superflow 
velocity is U = (0.07, 0) in dimensionless units. 



IV. RESULTS 



B. Magnus force on the lattice 

An immediate consequence of the above expressions 
is that changing the sign of vorticity, i.e., the signs of 
the coefficients a u and cr v , changes the motion of a vor- 
tex (which now becomes an anti-vortex) in such a way 
that the projection of the vortex velocity on the direc- 
tion of U remains the same, while the orthogonal pro- 
jection changes sign. In other words, a vortex and an 
anti-vortex move symmetrically with respect to the su- 
perflow. In general, for coarse lattices and superflow that 
is not parallel to one of the main periods of the lattice, we 
have found that the motion does not have that property. 
We interpret this as a result of anisotropy introduced by 
the lattice and by the superflow direction. To account for 
anisotropy, the net force in Eq. (fT2*|l needs to be replaced 
by 

F = -MV -<t v zxV + <t u zxU , (18) 

where M is a symmetric matrix that can depend on the 
direction of U. Such a matrix has three independent 
elements, which now replace the single drag coefficient of 
the isotropic model. 

On the other hand, if the superflow velocity is along 
one of the main periods of the lattice, the vortex and 
anti- vortex do move symmetrically with respect to it. In 
this case, we can introduce the effective Magnus force 
and drag coefficients that are defined by Eqs. (|16fl and 
Q17[l. In what follows, the superflow is always oriented 
along the x direction, U = (U x ,0), and we present two 
types of results: one type is the ratios V x /U x and V y /U x 
themselves, which are directly measurable quantities, and 
the other is the effective force coefficients computed from 
Eqs. (T5J and JTJJ. 



Simulations with different total lattice lengths have 
been carried out, with similar outcomes. The ratio of 
the vortex velocity to the superflow velocity for square 
lattices of the same length 600, (in the rescaled length 
units), constant superflow velocity U = (0.07,0), and 
different lattice spacings is shown in Fig. Because the 
vortex has to overcome pinning in the lattice cells, it 
moves by detectable jumps. The data points were ob- 
tained by averaging the measured vortex velocities over 
long time intervals that begin some time after the start 
of the simulation. 

We see that when the lattice spacing is close to 1, in 
our dimensionless units, we obtain results expected for 
the GI limit: V x /U w 1 and V y /U w 0, i.e., VmU, which 
means that vortices go with the flow. For large spacings, 
i.e., in the discrete limit, the behavior changes radically. 
Now, V x /U w and V y /U ^> 1, i.e., vortices move per- 
pendicular to the current. Between the two limits, there 
is an interesting regime when V x — V y , corresponding to 
motion with a Hall angle of 45°. 

Let us see how these results are reflected in the pa- 
rameters of our phenomenological model. As discussed 
in subsect. IIII Bl we determine the effective force coeffi- 
cients using Eqs. (TSJ, i|17|) . The result is shown in Fig. [5] 
In the quasi-continuous limit, we obtain a v j a u « 1 and 
T)/ a u *C 1, meaning that the Magnus force coefficient is 
nearly equal to the Lorentz force coefficient, while the 
drag coefficient is very small in comparison. 

In contrast, in the discrete limit, the Magnus force is 
nearly absent. This is our main result: we have observed 
a gradual vanishing of the Magnus force as we go from the 
continuous limit to the discrete one. We note that, while 
the Magnus force changes gradually, the drag coefficient 
has a rather sharp peak at the beginning of the region 
where a v ~ r/, i.e., when the vortex begins to move at 
45° to the flow. 
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Grid Spacing 

FIG. 5: (Color online) The ratios of the drag and Magnus 
force coefficients to the Lorentz force coefficient, as functions 
of the lattice spacing. The length of the lattice is 600; the 
superflow velocity is U = (0.07, 0) in dimensionless units. 
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FIG. 6: (Color online) Power spectra of ip in the middle of the 
simulation for different values of the lattice spacing dX. The 
length of the lattice is 600, and the lattice size (the number of 
lattice points) takes values 64 2 , 128 2 , and 256 2 . The straight 
line is a k~ 4 power law. 



The simplest intuitive picture that might account for 
the vanishing of the Magnus force is as follows. As the 
lattice spacing increases, the depletion of the density in 
the region of large phase gradients [the 8n(V9) 2 coupling 
discussed after Eq. becomes smaller, so that the 

system approaches the particle- hole-symmetric limit, in 
which the Magnus force is absent »iP A small remaining 
force, caused by deviations from this ideal limit, can be 
overcome by a force exerted by the lattice, resulting in 
zero net Magnus force. 

To make sure that the variations of the observables 
with the lattice spacing are not due to variation in 
the population of short-wavelength modes, we check the 
power spectra. Fig. shows the power spectrum of the 
field ip in the middle of the simulation for three different 



values of the lattice spacing. From the plot, we infer that 
there are no major differences in the power spectra. 

It is interesting to further explore the surprisingly 
broad intermediate range of lattice spacing where, as seen 

1.4 
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FIG. 7: (Color online) Hall effect (V y /V x = tan#„) as a func- 
tion of lattice spacing dY for lattices with equal and unequal 
spacings in the x and y directions. 



in Fig. the Hall angle is close to 45°. We have found 
that this preferred direction is related to the geometry 
of the lattice, i.e., it remains along the diagonal of the 
unit cell even after we go from the square unit cell to a 
rectangular one (see Fig.0). 



V. CONCLUSION 

The main result of the present paper is that in the 
classical limit vortices in superfluids on coarse (but 
uniform) lattices, in the presence of a superflow parallel 
to one of the main periods, move perpendicular to the 
superflow. We interpret this result as a transition from 
the full Magnus force in the Galilean-invariant limit to 
vanishing effective Magnus force in the discrete limit, in 
agreement with the observed smallness of the Magnus 
force in JJAi£ Our results are based on direct numeri- 
cal simulations of the discrete superfluid (J2J and do not 
assume a priori any symmetry that might prohibit the 
Magnus force in the discrete limit. 

Another potential application of our results is vortex 
motion in cold atomic gases confined in optical lattices. 
If a sufficient degree of experimental control over param- 
eters in either of these systems can be reached, interme- 
diate points in the transition from the Galilean-invariant 
limit to the discrete limit may become observable. No- 
tably, in our simulations these intermediate points in- 
clude a somewhat counterintuitive regime where, for a 
broad range of lattice spacings, the average vortex veloc- 
ity is along the diagonal of the unit cell. 
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